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Abstract 

It has recently been realized that the Balitsky-JIMWLK equations have 
serious shortcomings as equations to be used in small- x evolution near the 
unitarity limit. A recent generalization of the Balitsky equations has been 
given which corrects these shortcomings. In this paper we present an equiv- 
alent discussion, but in terms of the JIMWLK equation where we show that 
a new (fourth order functional derivative) term should be included. We also 
present a stochastic version of the new equation which, however, has some 
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1 Introduction 



The problem of high energy (small x) evolution near the unitarity limit [1] is one of 
the most important and widely studied problems in high energy QCD. A good under- 
standing of this problem is necessary for understanding the energy dependence of the 
saturation momentum and other properties of the color glass condensate (CGC) [2]. 
This in turn is important for good understanding of such diverse quantities as the 
initial energy density produced just after a high energy heavy ion collision and the 
Q 2 and x-dependence of structure functions in deep inelastic scattering at small x 
and at moderate values of Q 2 . 

The main equations that have been used to deal with high energy scattering near 
the unitarity limit are the Balitsky-JIMWLK equations [3-6] and the Kovchegov 
equation [7]. The Balitsky equations are an infinite hierarchy of coupled equations 
expressing the energy dependence of the scattering of high energy quarks and glu- 
ons (represented by Wilson lines in the fundamental and adjoint representations, 
respectively) on a target. The Kovchegov equation results when the large- N c limit 
of the Balitsky hierarchy is taken along with dropping all correlations in the target. 
The Kovchegov equation is a (nonlinear) equation for a single function and, while 
no exact solution is known, the general properties of the solutions of this equation 
near the unitarity boundary are understood [8-22]. The JIMWLK equation is a 
functional Fokker-Planck equation which expresses small-x evolution of the target 
in a way which is exactly equivalent to the evolution of the scattering amplitudes 
given by the Balitsky equations. The JIMWLK equation is the basic equation for 
the CGC. 

Geometric scaling [23] is a general property of solutions to the Kovchegov equa- 
tion, and the exact (asymptotic) energy dependence of the saturation momen- 
tum [16-18] is also known for the Kovchegov equation. However, recently two of the 
present authors argued that while guaranteeing that overall scattering amplitudes 
satisfy unitarity limits the Kovchegov equation allows significant unitarity violation 
in intermediate stages of small-x evolution [24]. These authors then proposed a sim- 
ple prescription for suppressing the unitarity violations (a simple boundary limiting 
small size rare dipole fluctuations) without, however, realizing the generality of their 
prescription. Later, in Ref. [25], it was realized that the prescription was identical 
to a procedure used to study reaction-diffusion systems in statistical physics [26, 27], 
and that the prescription has a general validity for a certain class of observables such 
as the energy dependence of the saturation momentum. Again using the relationship 
between small-x evolution and the time dependence of reaction-diffusion processes 
in statistical physics an asymptotic scaling law for the scattering amplitude, strongly 
violating geometric scaling was derived. This scaling differs from, but is related to, 
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that found in Ref. [24] where fluctuations were not included in the purely mean field 
approach followed there. 

Thus, the Kovchegov equation misses some essential aspects of high energy evo- 
lution near the unitarity boundary, although it may still be accurate for large nuclei 
or at intermediate energies. In this light the results of numerical simulations [28] 
comparing evolution using the Kovchegov equation with that using the Balitsky- 
JIMWLK equations come as a big surprise, because no essential differences were 
found. While this is partly due to the initial (Gaussian) distribution used in Ref. [28] 
the deeper understanding came from Iancu and Triantafyllopoulos [29] who noted 
that the Balitsky-JIMWLK equation also misses essential aspects of unitarity con- 
straints. In the language of "Pomerons" , Balitsky-JIMWLK has Pomeron mergings 
but no Pomeron splittings. The authors of Ref. [29] then propose a generalization of 
the Balitsky equations which includes a stochastic part, thus allowing for Pomeron 
splittings, in addition to the usual Balitsky hierarchy. (In the context of Pomeron 
interactions this term, and its relationship to the dipole model, has been known for 
some time [31].) It is likely that the Balitsky hierarchy augmented by this new term 
contains the essential ingredients for small- x evolution near the saturation boundary. 

The discussion of Ref. [29] is given completely in terms of the Balitsky approach 
to high energy scattering. In this paper we give the corresponding (physically equiv- 
alent) discussion in the context of the JIMWLK approach. That is, we give a new 
equation, which includes the original contribution due to JIMWLK along with an 
extra term, the resulting equation having the form of a functional differential equa- 
tion for the CGC weight function W[ct\. The JIMWLK equation, given below in 
Eq. (2.1), is a functional Fokker-Planck equation which has up to two derivatives 
of W[a] with respect to a but, in addition, has multiplicative powers of a such 
that there are always at least as many powers of a multliplying W[a] as there are 
a-derivatives acting on W[a]. The term having equal numbers of en's as derivatives 
of a generates the kernel for BFKL evolution while terms having extra powers of a 
correspond to "Pomeron mergings". The absence of terms having more derivatives 
than factors of a is the absence of "Pomeron splittings." The added new term has 
precisely four derivatives of a along with two factors of a. The exact nature of this 
term is determined by requiring that it correspond to one parent dipole splitting into 
two daughter dipoles in the dipole model. If the parent dipole is part of a BFKL 
evolution, and if the two daughter dipoles undergo BFKL evolution as generated 
by the original terms in the JIMWLK equation, then our new term can be viewed 
as a single Pomeron splitting into two Pomerons. In Section 2 we carry out the 
construction of this splitting term with the resulting modified JIMWLK equation 
given in Eq. (2.24). 
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Our construction of this new term is carried out in the weak field limit where the 
dipole model [35, 37] is valid. We have no argument to the effect that it correctly 
represents Pomeron splittings in a strong field environment. On the other hand this 
new term is only important, as compared with the original terms in the JIMWLK 
equation, in the weak field limit so that we believe (2.24) represents the essential 
physics of BFKL evolution and unitarity constraints. That is, we believe that (2.24) 
is a good "effective" equation for small- £ evolution. If one wished to include higher 
order BFKL-kernel effects in the JIMWLK equation one could expect corrections to 
our new term also to be necessary, but in the context of evolution using the lowest 
order BFKL kernel, (2.24) should be sufficiently accurate. Iancu and Triantafyl- 
lopoulos [36] have observed that (2.24) is exactly equivalent to their new equation, 
and we also expect that it is equivalent to the procedure of unitarization developed 
some time ago by Salam [30], at least in the energy regime where that procedure is 
valid. 

One of the very nice features of the JIMWLK equation is that it is possible to 
view it as a Langevin equation and hence use it as a basis for numerical simula- 
tions. Indeed, successful simulations have been performed by Rummukainen and 
Weigert [28]. This motivates trying to cast (2.24) in the form of a stochastic dif- 
ferential equation. In Section 3 this is done with, however, a few unusual aspects. 
In general a fourth order derivative will require a non-Gaussian noise, and this is 
illustrated by a simple example in Section 3.1. In Eq. (3.32) we give the rule for the 
evolution of a which, in a stochastic sense, is equivalent to (2.24). Eq. (3.32) involves 
three noise terms: ((z) is a fourth-order noise term, u a (x,uj) is a Gaussian noise 
term and £(x,u>) is a non-diagonal noise term whose correlator is given in (3.42). 
While we have been successful in casting our equation in the form of a stochastic dif- 
ferential equation the unusual form of the noise terms makes it unclear whether this 
form will be at all useful for numerical simulation. We do note, however, that the 
fourth-order noise terms have been discussed in the mathematical literature [33, 34], 
but we have not been able to find discussions of terms like our non-diagonal noise, 
and we are not certain that this term makes good mathematical sense. 



2 Weak field extension of the JIMWLK equation 

Iancu and Triantafyllopoulos [29] have recently recognized that the JIMWLK equa- 
tion [4-6] - the basic equation of the Color Glass Condensate formalism - includes 
Pomeron mergings but not also Pomeron splittings which are important at low gluon 
densities (weak fields) in the wavefunction of a target. In this Section we extend the 
JIMWLK equation in the weak field regime by adding the Pomeron splittings. 
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2.1 JIMWLK equation 

The JIMWLK equation [4-6] reads 

w w ^=lL sd®^ v) id<vi WrM - (21) 



where Y = Inl/x is the rapidity of the small- x gluons, J = J d 2 x (x, y denote 
transverse coordinates), a a is the color gauge field radiated by the color sources 
in the target, and the functional WVfa] is the weight function for a given field 
configuration. The kernel rj ab (x, y) is a non-linear functional of a, 

r] ab (x, y) = - I K(x, y, z) (1 - V^) - VjV y )* , (2.2) 

as it depends on the Wilson lines in the adjoint representation V and V\ 

Vi[a] = Pexp^ J dx-a a (x-,x)T a ^j , (2.3) 

where P denotes the path-ordering along the light-cone coordinate x~ and 

^ x >y> z >- (2tt)2 [x-z)\z-yY ■ [ZA) 

The average of a generic operator O is defined as 

(0) Y = J D[a] W Y [a] 0[a] (2.5) 
and its evolution with respect to Y, 

follows by first using (2.1) and then integrating twice by parts in (2.5). 

As an example, consider the scattering of a single dipole off a target. In the 
eikonal approximation, the T-matrix is given by 

T(x, y ) = 1-1- tr(v\x)V(y)) , (2.7) 

where V\x) (V(y)) is a Wilson line in the fundamental representation (T a — > t a in 
(2.3)) which describes the scattering of the quark (antiquark) off the color field in 
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the target, and the color trace divided by N c is the average over color. Inserting the 
above expression in Eq. (2.6), one obtains 



d /rr( a s N c f (x-y 

w (T(x,y)) Y = ^rj 



|2 



z) 2 (z - y) 2 



- T(x, y) + T(x, z) + T(z, y) - T(x, z)T(z, y)j y . (2.8) 

In this non-linear equation, the scattering of a single dipole {T(x, y))y on the l.h.s. 
is coupled to the scattering of two dipoles off the target (T^) = (T(x, z)T{z, y))y 
on the r.h.s., and, in general, the equation for (T^) will involve also (T < - Ar+1 - ) ). 
So, Eq. (2.1) generates an infinite hierarchy of coupled non-linear equations. The 
same set of coupled equations has also been derived by Balitsky [3] in a frame where 
only the projectile evolves with increasing rapidity. Taking the large N c limit of the 
JIMWLK-Balitsky hierarchy (non-dipolar terms become neglegible) and dropping 
the correlations in the target ({T(x,z)T(z,y)) — > (T(x,z)) (T(z,y)}), one obtains 
the Kovchegov equation [7] which is a closed equation for the T-matrix and, thus, 
much easier to handle than the coupled JIMWLK-Balitsky equations. Although no 
exact solution is known, the general properties of the solution to the Kovchegov 
equation are well understood [8-19]. 



2.2 JIMWLK equation in the weak field limit 

In the weak field limit (a ~ g< 1), i.e., for low gluon density in the target, one can 
expand the kernel rf h in Eq. (2.2) in powers of g, 



f }C(x,y,z)[a c (x)-a c (z)}[a d (y)-a d (z)} (2.9) 

J z 

which when inserted in (2.6) yields 

(2.10) 



0{O)y g 2 (T<T d ) ab 



dY 2tt 



where 



For the T-matrix (2.7) in the leading power in g, 

To(x,y) = ^- (a a {x)-a\y)f + 0{g\ (2.12) 



equation (2.10) leads to the dipole version of BFKL equation [35] 



(2.13) 

Note that in the weak field limit the JIMWLK equation (2.10) conserves the number 
of fields: two a's are created by a 2 and two are destroyed by (5 /5a) 2 . The BFKL 
equation (2.13) is an example which shows that the evolution of a two-point cor- 
relation function (l.h.s) involves only two-point correlation functions (r.h.s). More 
generally, in the weak field limit, the evolution equation for an N-point correlation 
function involves only N-point correlation functions. In the language of perturbative 
QCD, the conservation of the number of fields corresponds to processes in which no 
Pomeron splitting or Pomeron merging in the t-channel is allowed. 

The JIMWLK equation describes gluon (Pomeron) mergings in the hadron (nu- 
cleus) that are very important at high gluon densities (strong fields). This happens 
through the non-linearities in the JIMWLK equation: an N-point correlation func- 
tion on the l.h.s. of (2.10) can involve in general all the M-point functions with 
M > N on the r.h.s. once the higher orders in a (a n with n > 2) are taken into 
account in the expansion of rj ab in (2.9). Gluon mergings reduce the number of fields 
in a correlation function. Fig. la shows the merging of 4 gluons into 2 gluons (or a 
4-point function reduces to a 2-point function). A further BFKL evolution turns the 
4 to 2 gluon merging into a 2 to 1 Pomeron merging via the triple pomeron vertex. 




(a) GO 



Figure 1: Target evolution: (a) gluon merging, (b) gluon splitting. 

It has recently been recognized [29] that the JIMWLK equation does not provide 
gluon (Pomeron) splittings that turn out to be important for low gluon densities or 
weak fields in the target. This is because the JIMWLK equation (2.1) does not 
include higher order derivatives (S n /Sa n ) with n > 2 which in the weak field limit 
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(where rf h oc a 2 ) would generate gluon splittings. For instance, the JIMWLK 
equation does not include a 2 (S 4 /Sa 4 ) which would generate a splitting of 2 gluons 
into 4 gluons (a 2-point correlation function rises to a 4-point correlation function) 
as shown in Fig. lb. Thus, for low gluon densities in the target, there are no gluon 
number fluctations (Pomeron splittings) in the JIMWLK equation. 



2.3 Extension of the JIMWLK equation in the weak field 
region based on the dipole model 

In the next steps we will extend the JIMWLK equation (2.1) in the weak field regime 
by adding the gluon (Pomeron) splittings. We include the gluon number fluctuations 
in the JIMWLK equation based on the dipole model [35, 37] which is designed to 
properly describe the dilute region of a hadron wavefunction. The dipole model, on 
the other hand, has so far failed to accommodate gluon mergings, important at high 
gluon densities, which are a specialty of the JIMWLK equation. 

To construct the weight function HV[a] of a target in the weak field limit based 
on the dipole model, we first need to know the weight function for a single dipole 
Wy [ot] ■ The latter can be obtained by requiring that it gives the well-known result 
for dipole-dipole scattering in the 2-gluon-exchange approximation when used to 
evaluate the averaging for the To-matrix (2.12) according to (2.5). The following 
weight function for a single dipole with coordinates (x , y ) [38] 

WrM = { i + ^ j v g(«K »)g(«l«.. ») fo „ (u )L«(„) } *M ■ < 214 > 

gives the desired result for dipole-dipole scattering, 

(T (x,y)) = J B[a] W Yo [a]T (x,y) 



9 4 N g 
8N 2 

8N 2 



g(x\x ,y ) - G(y\x ,y )j 

J_ ln (x-y ) 2 (y-x ) 2 ^ 
4tt (x - x ) 2 (y - y o y 



(2.15) 



with the definitions 



X 2 jl 2 



/d k e 1 
(2^)2 = 4^ ln 

(2.16) 

and N g = N 2 — 1 . The infrared cutoff /i allows us to write down the propagator but 
it is harmless as it cancels in the difference in Q. 
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The weight function for N dipoles in the dipole model (DM) reads 



•*H°1 = E/*V Wn{H^.MM ^ } «[«] ■ (2.17) 
where 

ft(ti)=e(ti|ii-i,Xi) (2.18) 

with £Cj_i and £Cj the positions of the quark and the antiquark in the ith dipole. The 
shorthand notation Pn(Y) denotes the probability to find a system of N dipoles 
with the transverse coordinates (xq, Xi), (xi, x 2 ), ...(xn-i, Vo) via the evolution 
of an original dipole with coordinates (xq, yo) up to rapidity Y. The phase space 
integration is dT N = d 2 Xi d 2 x 2 ...d 2 x N _i and A = g 2 /(2N c ). 

The evolution of the weight function (2.17) with increasing rapidity Y is given 
by the following linear evolution equation 

QY Y M 2n 2 J x ^ z (x-z) 2 (z-y) 2 

a r 5 2 

~ + 2 L G{U ^ V)G{V ^ V) 5<*{u)8<*{v) ) 
A f 

A f 

x(l + - / g(u 2 \z,y)g(v 2 \z,y)—— — ) 
X (^N^> VlVla c (x)a c (y)Wr[*} ■ (2.19) 



This equation is a sum of two terms: a virtual term (loss term) and a real term (gain 
term) as we know them from the dipole model. The term V| a(x) a(x) plays the 
role of an annihilation operator for a dipole while (1 + A/2 J Q 2 (d / da) 2 ) corresponds 
to a creation operator for a dipole. In the virtual term a dipole is destroyed and its 
coordinates are transfered to the new one which is created. In the real term a dipole 
is converted into two new ones with the right coordinates and splitting factors. 

To explicitly show the annihilation and creation of dipoles, consider the evalua- 
tion of (2.19) with an arbitrary function O[o\ 

Z S /D W («P!W (2,0) 
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J will be zero unless one takes the functional derivatives from one of the dipoles in 
(2.17) and integrates them by parts in (2.20) on the ct(x) a(y) term. Using 



A 
2 



Qi{u)Gi{v) 



V 2 x V 2 y a c (x)a c (y) 



5a a (u)5a a (v) 

= XN g [S(x - Xi-x) - 5{x - Xi)\ [5{y - x^) - 5(y - Xi)\ 



one gets 

d_ 

dY 



AT N 

n,,.\ c ^ 

i=l 



2tt 



[Xi Xi — \ 



X 



+ 



X 



i_i,Xj,z ( x i—l z) 2 (z Xi) 2 

A f 6 2 
Qiu^xt-!, z)G(v 1 \x i _ 1 , z) 



(2.21) 



U1,V1 



(l + -/ </(ii2|s,a:i)0(«2|z,a:i 



5a a (tti)5a a ( , Ui)- 
5 2 



' U2,V2 



5a b (u 2 )5a b (v 2 ) 
5 2 



) 



oo „ AT - „ „2 

xE/^ wn{Hi/»w M , w ,) ) <h 

where dT' N means that the a3i_i and Xi coordinate integrations are missing. By 
comparing this equation with (2.19), we now explicitly see that the evolution is the 
same one known from the dipole model: When increasing the rapidity in one step, 
the ./V-dipole system either survives as it is (virtual term) or it evolves by splitting 
one of the dipoles into two new dipoles (real term). In the latter, the kernel (times 
a s N c / (2ir 2 )) in the first line in (2.21) is the right factor to turn P N into Pjv+i while 
the product of the terms in the third and fourth line are what is needed for converting 
the original ith dipole into two new dipoles (cCj_i, z) and (z, aj»). 

One can rewrite the evolution equation (2.21) also in the form 



dW$ M [a\ 
dY 



N=l 



dTN 



d P N {Y) 
dY 



N \ r s 2 



(2.22) 

which is consistent with Eq. (2.17). The identification of (2.22) with (2.21) gives 
the evolution equation for the probabilities (see [38]) 



dPN^Yjxi, ...,a;jv_i) a s N c 
dY ~ ' 2ir 2 



N 



i=i 



kCj_i - x^ 



(aui_i - z) 2 (xi 



P N (Y\x 1 , ...,x N -i) 
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. ®sN c [Xi-i ~ X i+1 ) 

where P/vC^ |cci, ccjv— 1) denotes the probability for iV-dipoles with coordinates 
(ceo, £Ci), (xi, X2), ...(xn-i, Vo) generated by the evolution of an original dipole with 
coordinates (x , y ). The first term is a loss term which describes the emission of one 
gluon from the original iV-dipole system while the second term is a gain term which 
describes the formation of the iV-dipole system via the splitting of one dipole in an 
original (N — l)-dipole system. It is the loss term which takes into account gluon 
(Pomeron) splittings missed in the JIMWLK equation [29,30]. Furthermore, it is 
easy to show that the probability is conserved by (2.23) [38]: J2f =1 J dT N P N (Y) = 1. 

Let's now compare the evolution equation (2.19) which follows from the dipole 
model with the JIMWLK equation in the weak field limit. The three QQ terms in 
(2.19) correspond to ordinary BFKL evolution while the QQQQ term allows Pomeron 
splitting. Since the ordinary BFKL terms were already properly included in the 
JIMWLK equation, only the final, QQQQ, term has to be added to the JIMWLK 
equation. Thus, the main result of this work reads 

N 2 f /A\ 2 (x-v) 2 

" S% y w (2 J (x - »)»("-„)» e(uAx - z)SMx - »>«*l*. ») 

A A A A \ 

x f 7 x f 7 x f r, s f r, x V' Vl a c (x) a c (y) W Y [a] . (2.24) 
5a a (wi) 5a a («i) 5a h (u 2 ) 8a b (v 2 ) y J 

Our derivation of the Pomeron splitting term, QQQQ term in (2.24), is done in the 
weak field limit. Therefore we cannot say that it also correctly describes Pomeron 
splittings in the presence of strong fields (a ~ l/g), or high gluon densities, in the 
target. However, since the Pomeron splitting term is important only in the weak field 
limit (a ~ g) as compared with the original JIMWLK term, we believe that (2.24) is 
a good "effective" equation for small-a; evolution which includes the essential physics 
of BFKL evolution and unitarity constraints. Our Pomeron splitting term should be 
accurate enough as long as the evolution involves only the lowest order BFKL-kernel. 
Corrections to the new term are expected if higher order BFKL-kernel effects are 
included in the original JIMWLK term. 

The extended JIMWLK equation (2.24) includes Pomeron mergings and Pomeron 
splittings and, thus, through iterations, Pomeron loops. This equation is equivalent 
to the new equation derived by Iancu and Triantafyllopoulos [36], and it should also 
be equivalent to the procedure of unitarization developed by Salam [30] in which 
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the gluon number fluctuations are properly treated. We expect that the solution 
to Eq. (2.24) will confirm the results of the recent works [24,25]. In these works it 
has been shown that due to the fluctuations the rate of growth of the saturation 
momentum (the "saturation exponent") considerably decreases and the geometric 
scaling [23] property is strongly violated as compared to the numerical results of the 
Balitsky-JIMWLK equation and the Kovchegov equation given in [28]. The latter 
equations miss the fluctuations or the Pomeron loops. 

Below we write the extended JIMWLK equation (2.24) in a form which is more 
appropriate for its stochastic interpretation given in the next Section [6,32]: 



Nl f f\\ 2 (x-yf 



G{ui\x, z)Q(vi\x, z)Q(u 2 \z, y)G(v 2 \z, y) 



8ir 3 N q Jw \2J (x-z) 2 (z-y) 2 

ui,ui,V2,va 

x f 7 \ F 7 x f t \ f t \ V* Vl a c (x) a c (y)] W Y [a] , (2.25) 
5a a (wi) <5a a (vi) 5a b (u 2 ) 5a b (v 2 ) x y v 1 yy> J 1 1 v 1 

with 

The original JIMWLK equation (first two terms in (2.25)) is the familiar functional 
Fokker-Planck equation [6, 32] which describes a diffusion process with a a (x) giving 
the drift term and rj a ' b (x,y) the diffusion coefficient. A Langevin equation for the 
fields a(x) can be written down which is equivalent to the original JIMWLK equa- 
tion [32]. The extended JIMLWK equation (2.25) is more complicated as it contains 
also the fourth order derivatives with respect to a (Pomeron splitting) in addition 
to the familiar Fokker-Planck equation. In the next Section we show how to obtain 
a stochastic differential equation also for the extended JIMWLK equation. The 
hope hereby is that one may use the stochastic differential equation for numerical 
simulations of (2.25) as done for the original JIMWLK equation in [28]. 



3 Interpretation of the evolution equation in terms 
of noise 

A functional equation like the JIMWLK equation although written in an elegant 
compact form (not so compact after adding our extension for a dilute color medium) 
is hard to use in actual computations. But this obstacle can be surmounted via the 
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observation that JIMWLK resembles very much a diffusion equation [32] which 
describes the time evolution of molecules undergoing Brownian motion. The latter 
has the interpretation of particles subjected to random acceleration in the medium, 
in other words a Langevin interpretation of the motion of the particles. In the 
case of JIMWLK, it is the light-cone gauge fields a a (x) that acquire a Gaussian 
random noise. In the previous section it was shown that in the case of a weak 
field a fourth order differential term was required to fully take into account all 
the important contributions at the leading log level. Unfortunately the extended 
JIMWLK equation ceases to be a functional diffusion equation of the Fokker-Planck 
type and the Langevin description is ruined. In this section we aim to restore the 
stochastic description even though the Langevin interpretation is mostly irreversibly 
lost. This requires a multi-random noise description involving both Gaussian white 
noise as well as other higher order noise. We will first show the basics of how a purely 
fourth order differential equation can arise out of a non-Gaussian random noise 
introduced at the fundamental level of the equation of motion. Then progressing 
in increasing complexity a Langevin term will be added which is well known to be 
responsible for the second order diffusion equation. Thus one will eventually have 
a multi-noise stochastic differential equation as the equation of motion from which 
the extended JIMWLK equation can be recovered. 

3.1 Illustration with some simple toy models 

To demonstrate the essentials of how to arrive at a fourth order differential equation 
using random noise, we consider a system of particles whose motions, similar to 
molecules undergoing Brownian motion, are Markov in nature. That means what 
happens to the system of particle next in time depends only on the present state 
and not in any way on the history of this system of particles. As we will see, 
this crucial condition simplifies the probabilistic description tremendously when the 
probability density distribution at one time can be related to a previous time in a 
relatively simple manner. Thus it is tempting to try to stay within the framework 
of Brownian motion, however we shall see that this is not really possible. 

3.1.1 A 1+1 dimensional system 

To start off, let's consider the probability distribution, P(x,t), of our system of 
particles which is a function of position x and time t. The motion of the particles 
is taken to obey the stochastic equation 
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where a(x) is the deterministic mean drift of the particles, 7(2;) is some function 
of x, and ((t) is a non-Gaussian random noise which has a normalized distribution 
(suppressing for the moment the t dependence) 

d ^NG(0 = / NG (C)dC (3.2) 

so that 

|/ NG (C)dC = l. (3.3) 

Using the notation 

( (O) NG = / ^(0/ NG (C)dC (3.4) 
to denote an average over the distribution / NG (C)> the random noise obeys 

(0 NG = (C 2 ) NG = (C 3 ) NG = , (3.5) 
and is only nonzero first in the four ( correlator 

(C 4 ) NG = c NG (3.6) 

where c NG is a number. This is the reason that ( is referred to as a fourth order noise 
in the literature [33, 34]. In that sense the usual Gaussian noise is second order noise. 
Note that the distribution / NG (C) cannot be written down simply as a function of 
some compact form but must be expressed in terms of a Fourier transform. The 
explicit expression / NG (C) serves little purpose in our discussion here and will not 
be written down. What is important here are the correlators given above, the first 
nonvanishing correlator is the one with four £'s whereas it is the two noise correlator 
in the case of the Gaussian noise that is the first nonvanishing one. 

Given that the particle motion is Markov in nature, the probability distribution 
of the particles P(x,t) must depend on P(x',t') where the particles currently at x 
arrived there by previously first reaching x' at an earlier time t' . Bayes' theorem 
states the conditional probability P(A\B) relating the probability P(A) of an event A 
happening given that another related event B has already happened with probability 
P(B) must satisfy 

P(A) = P{A\B)P{B) . (3.7) 



It follows that 



P(x,t) = J dx'P(x,t\x',t')P(x',t') . (3i 



For a small time increment St so that t' = t — St, the particle motion, from Eq. (3.1), 
changes by 

Sx(t) = x - x' = \a(x') + (7(x')) 1/4 C(*)l St ■ ( 3 - 9 ) 
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Since ( is a stochastic variable given the present position x, the previous x' is 
therefore non-deterministic but has a distribution. The conditional probability 
P(x,t\x',t') is therefore given by 



P(x,t\x',t') = J dP NG (C(t))5(x- a ;'-[ ( 7(x') + (7(^)) 1/4 CW]^). (3. 



10) 



In the limit that 5t becomes very small, one can Taylor expand the 5-function in 
powers of 5t, but one has to keep terms up to the fourth power of St if they are 
accompanied by a corresponding power of ( 



P(x,t\x',t') 



dP NG (CW) - x') - S'(x - x') [a(x>) + ( 7 (a;')) 1/4 C(0 



5t 



+ l -5'\x-x>)[a{x') + { 1 {x'))^at) 2 

~5"'(x-x') [a( a ;') + (7^ , )) 1/4 CW 
+ l5'"\x-x')[a{x<) + { 1 {x')yl\{t) 



5t 2 
5t 4 



+ 



...}. 



(3.11) 



Here we have written d5(x — x')/dx = 5'{x — x') = —d5(x — x')/dx'. In view of 
Eq. (3.5) and Eq. (3.6), the above expansion can be simplified to 

P(x, t\x', t') ~ 8{x - x') - 5'{x - x') <j{x') St + ^c NG 5""(x - x') 7 (x / ) <5t 4 + . . . (3.12) 
Substituting this into Eq. (3.8) and expanding 



P(x', t-6t)~ P(x', t) - dP< ^^ 8t +... 

at 



(3.13) 



(3.14) 



as well, after some rearrangements and using c NG = — I /St 3 , one finally has 

dP(x,t) dcr(x)P(x,t) 1 d 4 -f(x)P(x,t) 
Ft ~ dx 4! dx 1 ' 

The negative sign of c NG 4 is necessary for the stable evolution of P(x,t). One can 
see that fourth order diffusion can be associated with a fourth order noise much in 
the same way that the Gaussian white noise gives rise to the second order diffusion 
equation. 

4 This c NG is given by the inverse of the third power of the small timcstcp St is by no means 

accidental. From Eq. (3.6) one can deduce that ( oc <5t~ 3 / 4 or dV = (dt oc dt 1 / 4 . The last would 
be a more or less standard expression for the non-Gaussian stochastic variable V in the language 
of stochastic calculus. Equally one can arrive at Eq. (3.12) by applying Ito's lemma for fourth 
order noise. As one can see the equation has been expanded up to order 5t 4 only in appearance. 
In reality one has only reached the linear order in St. 
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3.1.2 A 2+1 dimensional system 



Let's now consider a system that has an evolution equation for its probability dis- 
tribution that bears more resemblance to the equation in Sec 2. Since the CGC 
evolution equation exists both in the functional space which is really infinite di- 
mensional as well as the color space of QCD, we extend our system to two spatial 
dimensions. This should suffice to demonstrate the main features of the noise inter- 
pretation of our equation. In addition JIMWLK is basically a diffusion equation so 
there must be a stochastic term for the white noise as well. We will let the motion 
of the particles in this system obey 

x\t) = o*(x) + <»(x) v j {t) + ^(xyf'GMuKt) ( b (t) . (3.15) 

o l (x) is the deterministic drift as before, which is followed by the term of the Gaus- 
sian white noise z/- 7 with the tensor coefficient function e tj (x). v l has the usual 
normalized Gaussian distribution 

dV G (u) = f G (u)d 2 u, (3.16) 

J / G (i/)d 2 i/ = l. (3.17) 

and expectation 

(0(")) G = J 0(u)f G (u)d 2 u. (3.18) 
The white noise z/(i) of course has all the usual properties 

H G = (3.19) 

and 

<!/*!/% =c G *«. (3.20) 

The last term of Eq. (3.15) is the fourth order noise term but now with more com- 
plicated functions as its coefficient as well as a set of secondary random Gaussian 
distributed variables v\. The last have the normalized distribution 



which gives 



dP G (P a ) = / G (P a )d 2 P a , (3.21) 



(^l) G =^5 ab . (3.22) 



The form of Eq. (3.15) is peculiar but is necessary to arrive at the extended 
evolution equation. Proceeding in a similar fashion to the previous subsection with 
a small time increment 

5x\t) = x l - x' 1 = \a\x') + e ij (x') v j (t) + ( 7 (a0) 1/4 G ab (x')Pi(t) { b (t)]6t . (3.23) 
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Bayes' theorem gives us 

P{x,t) = J d 2 x'P{x,t\x't')P{x',t') , (3.24) 
with the conditional probability 
P(x,t\x't') = f dV G (v) ]JdP G (9 a ) dV NG (() 

J a 

x< j2 (a . _ x >_ [a i^ + e ij {x/) v i {t) + ( 7 ( a '))V*G afc (a / )P*(0 a(oi^) 

(3.25) 

which once again can be expanded or one can equally apply Ito's lemma. Bearing 
in mind all the vanishing correlators of the different noise, we have 

p( X ,t\x't') = i dv G (u) n d ^ G (^) d ^N G (c) 

J a 

{5 2 (x - x') - a'ix') d i 5 2 (x - x')5t 

+^e im (x')e jn (x')is m (t)is n (t)d i d i 6 2 {x - x') 5t 2 

+^l{x') G ab {x')G cd {x')Ge f {x')G gh {x') 9i(t)9i(t)P k e (t)9 l g (t) 

x(b(t)Ut)(f(t)( h (t) d^d k d l 5 2 (x - x') 5t 4 
+ ...}. (3.26) 

After making use of Eq. (3.20), Eq. (3.22), the modified version of Eq. (3.6), 

(( b (t)Ut)( f (t)( h (t)) = c NG 5 bd 5 df 5 fh , (3.27) 

the notation from [32], 

rj ij (x) = e ik (x)e jk (x) , (3.28) 

and not forgetting that the correlator of four Gaussian noise isn't zero but can be 
separated into a sum of products of two two Gaussian noise correlators, 

(m^ l 9 ) G = {m) G {»y g ) G + fe%<^> G + <^> G fe% , (3.29) 

one gets 

dP(x,t) _ da\x)P(x,t) | l d 2 ^(x)P(x,t) l d A A^ kl (x)P(x,t) 
dt dx % 2! dx l dxi 4! dx l dxWx k dx l ' 
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where c Q = 1/St has been taken 5 and 

A ijk \x) = 7 ( a! )G a6 (a;)G c6 (x)Ge 6 (x)G ff6 (x) 

x (5 ac 6 eg 5 ij 6 kl + 5 ae 5 cg 5 tk V l + 5 ag 5 ce 5 u 5^ . (3.31) 

This equation has all the basic features of Eq. (2.25) but some fine tuning is still 
required. This is the subject of the next subsection. 

3.2 Noise interpretation of the extended JIMWLK equation 

Back to our present problem, the JIMWLK equation is a functional equation of the 
integrated light-cone gauge field a a and the evolution is in rapidity Y rather than 
in time so one has to replace t by Y everywhere. In analogy to the equation for x % 
before, we write the stochastic equation 

JL a a (u) = a a (u)+ [ e ab ' i (u, Z y' i (z)+ [ p(u, x, w, z)P a (x, w)((z) yftfcri) 

J z J x,w,z 

(3.32) 

where borrowing again from the notation of [32] 6 the first stochastic term with the 
Gaussian noise v h > % has the tensor coefficient function 

e^(u, z) = -=) -^(1 - VlV z r (3.33) 

V47T 3 \u - zy 

and the product 

rf{u, v) = J e ac ' l (u, z)e bc ^(v, z) (3.34) 
from which the first deterministic term cr a (u) can be expressed as 

The last multi-noise stochastic term on the right hand side of Eq. (3.32) has a second 
but different Gaussian noise is a (x,w), the fourth order noise ((z) and a new noise 
£(jc, w ) which we call a non-diagonal noise for reason which will become clear later 
with a normalized distribution 



dP ND (0 = l. (3.36) 



5 In this case with = 1/St the white noise v l is of the order of v % ~ \j\JHt. In terms of 

the conventional stochastic variable that would be dW l — v l dt ~ Vdi as required from stochastic 
calculus. 

6 Note the slightly different notation of the a a of [32] which is equivalent to our da a /dY. 
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The coefficient function of this last term is 



Here we consider an evolution of the probability density Wy [oi\ of the color fields 
of the CGC by a small increment from Y' = Y — 5Y to Y 

W Y [a}= J V[a']Py\y,[a\a']W Y ,[a'] (3.38) 

during which the fields a a make the change 

Sa a = a a - a' a = 5Y ^-a a . (3.39) 

oY 

The conditional probability here is 

Py\Y>WW'\ = /n[ D ^G(^)]ri[ D ^G(^)][ D ^NG(C)][D^ND(0] 

a,i b 

x 5 N " (a -a'- 5Y da/dY) . (3.40) 

This is written slightly differently than before because of the functional aspect of 
Eq. (3.38). The measures of the noise [DP] here are now path integrals for the same 
reason, and da/dY is given by Eq. (3.32). Either using Ito's lemma or expanding 
the delta function in terms of SY, Wy[a'] in terms of 5Y as before and substituting 
in Eq. (3.38), the original JIMWLK equation is straightforwardly recovered by fol- 
lowing the steps in the previous subsection [32]. The term that is less trivial is the 
four derivative Pomeron splitting term which comes essentially from the last term of 
Eq. (3.32) raised to the fourth power or the collection of all the apparent 5Y 4 terms 
in the expansion of the delta function in Eq. (3.40). At this order of the expansion, 
each functional derivative 8/5a a is accompanied by ((z) whose correlator ensures 
that there is only one z coordinate via 

(C(z)C(z0C(^)C(O) NG = c NG 5\z - z')b\z' - z")5\z" - z">) , (3.41) 

where here c NG = — 1/SY 3 . This is the case as seen in Eq. (2.25), the new term 
has four Q functions sharing the same z as one endpoint to a pair of dipoles. This 
is the result of a one-step evolution from a single parent dipole. Each of the four 
Q(u\x,zYs coming from the fourth order expansion of Eq. (3.40) is a contribution 
to the field a(u) from a dipole with endpoints x and z. Or equivalently each dipole 
interacts with the target through two gluon exchanges, each of which is represented 
by a Q(u\x, z) therefore the four Q(u\x, z)'s must come in two pairs with each pair 
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sharing the same endpoints. The presence of the second Gaussian noise v a (x,w) 
correlator (u a (xi, w 1 )is b (x2, w 2 )i /C (x 3 , w 3 )v d (x 4i) 1/74) ) G , which can be rewritten in 
a similar fashion to Eq. (3.29), guarantees that this is the case. These mostly 
ensure that one inevitably arrives at Eq. (2.25) except the numerator of the BFKL 
kernel, which is (x — y) 2 . Naively one is required to put the quartic root of this 
expression in the stochastic equation for a a in Eq. (3.32) which is impossible since 
only one of x and y can appear in the expression. The other comes into Eq. (2.25) 
via the 5Y expansion and the v correlator discussed above. Thus to ensure that 
the appearance of the correct BFKL kernel, we gave the factor {{x — w) 2 } 1 ^ 4 to 
Eq. (3.32), introduced the very unconventional non-diagonal noise y/£(x, w) and 
required the latter to satisfy the correlation 

{£{x, w)£(y, w')) ND =5 2 (x- w')5 2 {y - w) . (3.42) 

One should note the unusual pairing of the arguments of the 5-functions. This 
is the reason for its nomenclature. Here two sets of coordinate pair have already 
been identified as x, w and y, w' thanks to the v a correlators so the square root has 
been removed. This correlation is somewhat counterintuitive because one encounters 
much more frequently the normal Gaussian noise, whose correlator is 

(£(*, w)t(y, w')) Q = 5 2 (x - y)5 2 (w - w') . (3.43) 

As it is the non-diagonal correlator requires that a squared quantity to be zero 

<£0r,<u,) 2 ) ND = O. (3.44) 

This demands a lot on the distribution itself. The only constraint on £(x, w) here is 
Eq. (3.36) and Eq. (3.42). Note that although ^/£(x, w) appears in the stochastic 
equation Eq. (3.32), one does not encounter correlators such as 

etc because of the presence of u a (x, w) and ({x). Such correlators are not only un- 
pleasant but also not very meaningful. We are unable to find such an unusual corre- 
lator in the mathematics literature but if the stochastic interpretation of Eq. (2.25) 
is to be preserved, one would have to require such a correlator to exist. Note that 
the vanishing of the average over a squared quantity similar to Eq. (3.44) is also 
required of the fourth order noise ( for which one knows that a distribution does 
exist. Providing one accepts that a distribution for £ might probably exist then 
with all the ingredients in place Eq. (2.25) can be reproduced from the multi-noise 
stochastic equation Eq. (3.32) by going through the steps as shown in this section. 
So once again although the extended JIMWLK equation for a dilute color dipole 
medium appears to be more complicated with a fourth order functional derivative 
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term so that the Langevin description is no longer an option, nevertheless one can 
preserve the stochastic interpretation through the introduction of higher order noise. 
The original goal for the Langevin description is to facilitate the computation of the 
evolution of high energy collisions with CGC through the JIMWLK equation. In a 
dilute medium the new term is as important as any terms in the original JIMWLK 
equation, however computationwise, it is not clear that the new stochastic equation 
here necessary for the new development would fulfill its original purpose since it 
might be difficult, if not impossible, to implement both the fourth order noise and 
the non-diagonal noise on a computer. Again provided a distribution exists, the last 
would likely to pose the most difficulties. Should this distribution turn out eventu- 
ally to be not possible then the stochastic description is really lost. We will leave 
that for future consideration and discussion. 



Acknowledgements 

We thank Edmond Iancu for many enlightening discussions on this general prob- 
lem and for explaining to us the ideas and details of his work with D. Triantafyl- 
lopoulos. We thank the latter for an elegant note showing the equivalence between 
their approach and ours. A. Sh. acknowledges financial support by the Deutsche 
Forschungsgemeinschaft under contract Sh 92/1-1. 



References 

[1] L. V. Gribov, E. M. Levin and M. G. Ryskin, Phys. Rept. 100 (1983) 1. 

[2] For reviews see: E. Iancu, A. Leonidov and L. McLarren in "QCD Perspectives 
on Hot and Dense Matter", J. -P. Blaizot and E. Iancu eds., Kluwer Academic 
Publishers (2002), hep-ph/0202270; A. Mueller in "QCD Perspectives on Hot 
and Dense Matter" , J.-P. Blaizot and E. Iancu eds., KluwerAcademic Publishers 
(2002), hep-ph/0111244; E. Iancu and R. Venugopalan, hep-ph/0303204. 

[3] I. Balitsky, Nucl. Phys. B 463 (1996) 99; Phys. Rev. Lett. 81 (1998) 2024; 
Phys. Lett. B 518 (2001) 235. 

[4] J. Jalilian-Marian, A. Kovner, A. Leonidov and H. Weigert, Nucl. Phys. B 504 
(1997) 415; Phys. Rev. D 59 (1999) 014014. 

[5] E. Iancu, A. Leonidov and L. D. McLerran, Phys. Lett. B 510 (2001) 133; Nucl. 
Phys. A 692 (2001) 583. 



21 



[6] H. Weigert, Nucl. Phys. A 703 (2002) 823. 

[7] Y. V. Kovchegov, Phys. Rev. D 60 (1999) 034008; Phys. Rev. D 61 (2000) 
074018. 

[8] E. Levin and K. Tuchin, Nucl. Phys. B 573 (2000) 833; Nucl. Phys. A 693 

(2001) 787. 

[9] N. Armesto and M. A. Braun, Eur. Phys. J. C 20 (2001) 517; Eur. Phys. J. C 
22 (2001) 351. 

[10] K. Golec-Biernat, L. Motyka and A. M. Stasto, Phys. Rev. D 65 (2002) 074037. 

[11] E. Levin and M. Lublinsky, Phys. Lett. B 521 (2001) 233; Eur. Phys. J. C 22 

(2002) 647; M. Lublinsky Eur. Phys. J. C 21 (2001) 513. 

[12] K. Golec-Biernat and A. M. Stasto, Nucl. Phys. B 668 (2003) 345. 

[13] J. L. Albacete, N. Armesto, A. Kovner, C. A. Salgado and U. A. Wiedemann, 
Phys. Rev. Lett. 92 (2004) 082001. 

[14] J. L. Albacete, N. Armesto, J. G. Milhano, C. A. Salgado and U. A. Wiedemann, 
"Numerical analysis of the Balitsky-Kovchegov equation with running coupling: 
Dependence of the saturation scale on nuclear size and rapidity," arXiv:hep- 
ph/0408216. 

[15] E. Iancu, K. Itakura and L. McLerran, Nucl. Phys. A 708 (2002) 327. 
[16] A. H. Mueller and D. N. Triantafyllopoulos, Nucl. Phys. B 640 (2002) 331. 
[17] D. N. Triantafyllopoulos, Nucl. Phys. B 648 (2003) 293. 

[18] S. Munier and R. Peschanski, Phys. Rev. Lett. 91 (2003) 232001; Phys. Rev. 
D 69 (2004) 034008. 

[19] S. Munier and R. Peschanski, Phys. Rev. D 70 (2004) 077503. 

[20] J. Bartels, V. S. Fadin and L. N. Lipatov, Nucl. Phys. B 698 (2004) 255. 

[21] G. Chachamis, M. Lublinsky and A. Sabio Vera, "Higher order effects in non 
linear evolution from a veto in rapidities," arXiv:hep-ph/0408333. 

[22] E. Levin and M. Lublinsky, "Balitsky's hierarchy from Mueller's dipole model 
and more about target correlations," arXiv:hep-ph/0411121. 

[23] A. M. Stasto, K. Golec-Biernat and J. Kwiecinski, Phys. Rev. Lett. 86 (2001) 
596. 



22 



[24] A. H. Mueller and A. I. Shoshi, Nucl. Phys. B 692 (2004) 175; "Small-x physics 
near the saturation regime," arXiv:hep-ph/0405205. 

[25] E. Iancu, A. H. Mueller and S. Munier, "Universal behavior of QCD amplitudes 
at high energy from general tools of statistical physics," arXiv:hep-ph/0410018. 

[26] E. Brunet and B. Derrida, Phys. Rev. E 86 (1997) 2596; Comp. Phys. Comm. 
121-122 (1999) 376; J. Stat. Phys. 103 (2001) 269. 

[27] For a recent review, see W. V. Saarloos, Phys. Rep. 386 (2003) 29. 

[28] K. Rummukainen and H. Weigert, Nucl. Phys. A 739 (2004) 183. 

[29] E. Iancu and D. N. Triantafyllopoulos, "A Langevin equation for high energy 
evolution with pomeron loops," arXiv:hep-ph/0411405. 

[30] G. P. Salam, Nucl. Phys. B 449 (1995) 589. 

[31] J. Bartels and M. Wusthoff, Z. Phys. C 66 (1995) 157; 

M. A. Braun and G. P. Vacca, Eur. Phys. J. C 6 (1999) 147; 

J. Bartels, L. N. Lipatov and G. P. Vacca, "Interactions of Reggeized gluons in 

the Moebius representation," arXiv:hep-ph/0404110. 

[32] J. P. Blaizot, E. Iancu and H. Weigert, Nucl. Phys. A 713 (2003) 441. 

[33] C. W. Gardiner, Handbook of Stochastic Methods, Springer Verlag, Second Edi- 
tion, 1985. 

[34] K. J. Hochberg, The Annals of Probability 1978, Vol. 6, No. 3, 433-458. 
[35] A. H. Mueller, Nucl. Phys. B 415 (1994) 373. 

[36] E. Iancu and D. N. Triantafyllopoulos, paper in preparation, private communi- 
cation. 

[37] A. H. Mueller and B. Patel, Nucl. Phys. B 425 (1994) 471. 
[38] E. Iancu and A. H. Mueller, Nucl. Phys. A 730 (2004) 460. 



23 



